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Abstract. We analyze the classical version of a plaquette orbital model that was recently 
introduced and studied numerically by S. Wenzel and W. Janke. In this model, edges of the 
square lattice are partitioned into x and z-types that alternate along both coordinate directions 
and thus arrange into a checkerboard pattern of x and z-plaquettes; classical 0(2)-spins are 
then coupled ferromagnetically via their first components over the x-edges and via their 
second components over the z-edges. We prove from first principles that, at sufficiently low 
temperatures, the model exhibits orientational long-range order (OLRO) in one of the two 
principal lattice directions. Magnetic order is precluded by the underlying symmetries. A 
similar set of results is inferred also for quantum systems with large spin although the spin- 1 ^ 
instance currently seems beyond the reach of rigorous methods. We point out that the Neel 
order in the plaquette energy distribution observed in numerical simulations is an artefact of 
the OLRO and a judicious choice of the plaquette energies. In particular, this order seems to 
disappear when the plaquette energies are adjusted to vanish at the ground-state level. We 
also discuss the specific role of the underlying symmetries in Wenzel and Janke's simulations 
and propose an enhanced method of numerical sampling that could in principle significantly 
increase the speed of convergence. 
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1. Introduction 

The physics of transition-metal compounds [20] and of the underlying effects such as 
colossal magnetoresistance [8] has recently spawned a number of spin models of high 
theoretical and practical interest [26, 36]. Among the common features of these models 
is that (1) the degrees of freedom are represented by quantum or classical spins residing 
at the vertices of a regular lattice, typically, the square or cubic lattice, (2) the interactions 
are nearest neighbor and ferromagnetic but (3) only certain components — or projections 
— of the spins are coupled over each edge of the lattice. The spin variables actually 
represent effective degrees of freedom (pseudospins); typically, occupation characteristics 
of a partially filled atom orbital. The interaction is also effective and it is arrived at by 
considerations of crystal fields mediated by interlaced atoms [25], or by appealing to Jahn- 
Teller distortions [22] (or both methods, via different routes, [39]). Systems of a similar 
nature have sprung up independently in the field of topological quantum computation, 
e.g., the Kitaev model [28, 27], and so besides practical incentives to develop a theory for 
the behavior of these models, there are also strong theoretical reasons to understand their 
possible technical implementations [15, 23]. 

The definition of the aforementioned class of models starts by partitioning all edges of the 
lattice into families indexed by some a; a generic edge in the a-th family is then denoted by 
(r, r')a- The Hamiltonian invariably takes the form 

H:=-EE J*sl a) sM (1.1) 

* {r,r')a 

(a.) 

with positive coupling constants, f a > 0. The explicit meaning of the projections S r is then 
a matter of what specific model one wishes to consider. 

Two examples of interest have been studied earlier: the orbital compass model (e.g., 
[24, 31, 13, 33, 34, 9]), where s£ , a. = 1, 2, 3, are the corresponding Cartesian components of 
the quantum spin and (r, r') a is an edge in the a-th lattice direction, and the 120-degree model 
(e.g., [32, 4, 5]), where the meaning of the edge (r, r') a is preserved but Sf , a = 1, 2, 3, now 
denotes the projections of the (three-component) spin S r onto the vectors 

i7 i: =(l,0), v 2 :={-\,%, -°s-= (1-2) 

(a) 

i.e., S r := v K • S r . Kitaev's model is defined similarly to the orbital compass model but the 
underlying graph is the honeycomb lattice. 

As usual, all models in the above class have a natural quantum version, where S r is a 
three-component spin operator — with a distinct irreducible matrix representation for each 
non-negative half integer — and a classical version, where S r is a vector a priori uniformly 
distributed on the unit sphere in R N (i.e., an O(N)-spin). 
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1.1. The plaquette orbital model 

Recently, an interesting variant of the orbital compass model has been proposed and studied 
by Wenzel and Janke [37]. Their model, which they termed the plaquette orbital model (POM), 
is most naturally defined over the square lattice Z 2 , although generalizations to higher 
dimensions are straightforward. The index ot takes only two values, a = 1, 2, and the spin 
projections are defined as follows: 

Si 1] := S x r and S? } := S z r , (1.3) 

where, in agreement with Wenzel and Janke's notation, (S*,Sj;) denote the Cartesian 
components of the vector-valued 0(2)-spin S r in the classical version while, in the quantum 
version, it denotes the corresponding pair of operators for the quantum spin. 



nxrxn 




nxrxn 

Figure 1. The configuration underlying the definition of the plaquette orbital 
model. Here the x-components of the spins are coupled over the red (solid) 
edges and the z-components are coupled over the blue (dashed) edges. 

What distinguishes this model from the previously discussed counterparts is the 
partitioning of the lattice: edges are designated as x and z-types in an alternating fashion 
along each line of sites (copy of Z) inside Z 2 so that all of the lattice faces (plaquettes) with 
the lower-left corner on the even sublattice of Z 2 contain only edges of one type. Thus one 
gets the Hamiltonian in accord with the general form in Eq. (1.1) with the edges of two types 
arranged into to a checkerboard pattern as in Fig. 1, 



U := ~h E S r S r' ~ h E S r S r» 

(rS)» (r,r') z 



(1.4) 
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with the x-components of the spins coupled over the x-edges and z-components over the 
z-edges. The signs of the coupling constants are immaterial as they can always be absorbed 
into a sign change of the corresponding component on one of the sublattices. 

In [37] Wenzel and Janke studied the POM numerically both in its classical and quantum 
version. For an order parameter they chose the plaquette energy 

C / C^l C«l I C«2C lt 2 I C K 3 C a 3 i C a 4 C a 4 \ /I c\ 

Of .— <^D r J r+ei -T D r >J r + e2 -T J r+ei J r+ei+e2 -T D r+ei D r+ei + e2 j, yi.D) 

indexed by the position of the lower-left corner of the plaquette, where e\ and ej_ are the unit 
vectors in the first and second coordinate direction and dl\, . . . , a.4 are either x or z depending 
on whether the corresponding edge is an x or z-type, respectively The simulations of 
[37] indicated a clear onset of Neel order in the plaquette energies at low temperatures 
for both the classical two-component (0(2) -spin) model, and the quantum, spin-l/2 model. 
Explicitly, in one such state the plaquette energies on the x-plaquettes are low but those 
on the z-plaquettes are generally high and in another state the roles of x and z-plaquettes 
are interchanged. (The mixed-type plaquettes end up with an intermediate energy in both 
circumstances.) However, no specific conclusion was attempted for other possible ordering 
mechanisms (magnetic order, orientational order, etc). 

1.2. Goals and outline 

In the present paper we wish to address the true nature of the phase transition in the POM 
by means of rigorous mathematical (i.e., analytic) methods that draw from earlier work 
on models of this kind [32, 4, 1, 5]. We will predominantly focus on the classical 0(2)- 
spin version of the model as the quantum spin- 1 /2 version poses technical difficulties that 
we do not yet know how to overcome. Nonetheless, thanks to the general theory [5], the 
conclusions for the classical system permit straightforward extensions to quantum systems 
once the magnitude of the quantum spin (i.e., the quantity S such that S r ■ S r = S(S + 1)) is 
sufficiently large compared to the inverse temperature squared. 

A key conceptual difference with the approach of Wenzel and Janke [37] is that instead 
of going via plaquette energies, we directly attempt (and succeed in) proving orientational 
long-range order (ORLO) of the spins in one of the two principal lattice directions. We then 
argue that the Neel ordering in plaquette energies found in [37] is actually an artefact of 
this ORLO: indeed, it is a direct consequence of the alignment of the spins along coordinate 
axes, the 2-periodicity of the interaction and the fact that the plaquette energies were not 
normalized to vanish in the ground states. In fact, with such normalization the Neel order 
seems to disappear altogether. 

The remainder of this paper is organized as follows: In Sect. 2 we will discuss the ground 
states of the classical Hamiltonian and then state our principal results concerning the ORLO 
in the pure and diluted systems. In Sect. 3-5 we provide full and reasonably self-contained 
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mathematical proofs of these results. In Sect. 6 we discuss connections with the numerical 
findings and propose a potentially more efficient way to obtain samples of equilibrium 
configurations in this model. 



2. Rigorous results 

2.1. Ground states 

Our discussion of the results opens up with the description of the ground states. We will 
focus on the situation in finite volumes with periodic boundary conditions. Specifically, let 
Tjv be the N x N torus — obtained by periodizing the square {0, ...,N — 1} x {0, . . . , N — 1} 
— and assume that N is even to reflect the natural period-2 nature of the interaction. Let Hn 
denote the Hamiltonian on T;v which we define by (1.4) with the edges (r, r') a restricted to 
nearest neighbor pairs (of the proper type) on Tjy. As is common, we call a configuration 
S = (S r ) r6 x N aground state of %n if 

H N (S) =mmH N {S'). (2.1) 
s' 

Here we note that the minimum is attained — and a ground state exists — by the sheer fact 
that Hn is a bounded and continuous function on a compact metric space. The issue is how 
many ground states there are and how they can be concisely described. 

A common feature of all models of the type (1.1) is abundance of symmetries with respect 
to simultaneous flips of (specific) components of the spins. An indisputable advantage of 
the POM over the other systems is that here the flips can be made locally. Explicitly, for r 
with both coordinates even, i.e., the site designating a lower left corner of an x-plaquette, 
let tp r (S) be the configuration defined by 

[q> (S)] x ■= { ~ SXf ' i° rr/ = r,r + e lf r + e 2 ,r + e 1 + e 2/ ^ 
I S*/ otherwise, 



and 



Wr(S)t := S z r , (2.3) 



For r with both coordinates odd (lower left corner of a z-plaquette), the map q> r is defined 
in a similar manner — with the z-components reflected instead of the x-components. A 
moment's thought now shows that 

n N (<p r (S)) = H N (S), (2.4) 

i.e., S i y q> r {S) is a symmetry of the model. 
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Theorem 2.1 Suppose ]\ = ]i> 0. Then every ground state of 'Hn can be obtained from a constant 
configuration, S r = e, for some unit vector e G R 2 , by successive applications of a subset of the 
maps ((p r ). 

tf h > h > ®' th en tdl ground states arise (via applications of (<p r )) from S r = e\, while if 
h > h > 0, then all ground states arise (again, via (cp r )) from S r = ej_. 

This statement is a precursor of the orientational LRO that we will establish for low 
(but non-zero) temperatures. The key problem there will be the fact that the ground-state 
degeneracy significantly increases at the symmetry point ]\ = ]i- This appears to be a 
common feature for all models covered by the equation (1.1). 

We also remark that in infinite volume (i.e., the model defined on all of Z 2 ), the structure 
of ground states is considerably more complicated. (In infinite volume, a ground state is 
a configuration whose any local change will result in a non-negative change of energy.) It 
is still true that any configuration obtained from constant configurations by means of the 
maps (q> r ) is a ground state. However, further ground states can be constructed by imposing 
linear (or other) interfaces. 



2.2. Orientational order 

We proceed to discuss our results for positive temperatures. Let jS := ^ denote the inverse 
temperature. In the canonical ensemble, the spins on are distributed according to the 
Gibbs measure Unr that is defined by 

e -pn N (s) 

u N/fi {dS) := — J] v(dS r ), (2.5) 

where v denotes the uniform (Haar) measure on the unit circle in R 2 normalized, for later 
convenience, to y/2n. The normalization constant Zjv,^ is the partition function. 

In order to formulate the existence of a long range order in a mathematically precise 
way, one often considers Gibbs measures directly in infinite volume — i.e., as measures on 
infinite configurations (Sj-) r6Z 2- These are defined by means of the so-called DLR-condition 
[11, 12, 30] stating that the conditional distribution in any finite A C Z 2 given (S r ) r6Z 2\ A 
takes the above form with Ti^ replaced by the Hamiltonian in A under the boundary 
conditions (S r ) r6Z 2\ A - 

A standard reference to the corresponding mathematical theory of these measures is 
Georgii [19]. We will now list the aspects that have a bearing on our problem. It is a 
standard fact that, for compactly-supported spins, infinite-volume Gibbs measures can be 
extracted (as weak limits) from sequences of finite-volume, or even torus, Gibbs measures. 
In particular, estimates on u^ ^ that hold uniformly in N readily yield corresponding 
estimates for infinite volume limits of {^n,^}- 
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In light of the period-2 nature of the interactions, a (Gibbs) measure ji will be called 
translation-invariant if u(c r a (A)) = }i(A) for any event A, where 

[o- a (S)] r :=S r -2e K , oc = 1,2, (2.6) 

are the shift operators. This measure is said to be ergodic if }i{A) is either zero or one for 
any event A that is invariant under the translations above, i.e., such that a~ 1 {A) = A for 
oc = 1,2. Ergodic measures have the distinguished feature that all block averages converge 
to the corresponding expectations (i.e., expected value represents the average value in 
sufficiently large boxes) and that a typical sample from these measures thus has a fairly 
homogeneous structure at sufficiently large spatial scales. 

We begin by ruling out magnetic ordering at any parameters of the model. 

Theorem 2.2 Let ]\, ]i > and /3 > 0. Then 

E f( (S r )=0 (2.7) 
for all infinite-volume Gibbs states u. 

The absence of magnetic order is of course a consequence of the symmetry S h- \ q> r (S). 
This result is prototypical for all models of the kind (1.1). It does not, however, exclude the 
existence of an orientational long-range order. 

Theorem 2.3 Suppose Ji = h > 0. Then for each 6 > 0, there is /3o = (5o(6) and a sequence (en) 
with £n — > such that for each /3 > /3o(^), 

mp( E [Sr] 2 >\^N\(l-S)) >l-e N , cc = x,z. (2.8) 

In particular, for each j6 > Po{$), there exist two translation-invariant, ergodic infinite-volume Gibbs 
measures, ji x ^ and u 7 ^, such that 

E ?( «([S?] 2 ) >l-5, oc = x,z. (2.9) 

In addition, for any > (5o(5) and any translation-invariant, ergodic infinite-volume Gibbs state u 
at inverse temperature f>, 

E F ([S?] 2 ) e [0,5]U [1-6,1], oc = x,z. (2.10) 

Notice that (2.8) states that typical configurations sampled from u^^ on a large torus 
have a majority of all spins aligned either along direction ±e\ or direction ±^2 (with both 
orientations equally present thanks to the symmetries (2.4)). The inequality (2.9) is a version 
of this fact in infinite volume and it manifestly demonstrates the occurrence of a phase 
transition. The statement (2.10) in turn implies that only the ground states oriented along 
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the principal axes are stable under thermal perturbations. The infinite degeneracy at the 
symmetry point for /3 = oo is thus reduced to a two-fold degeneracyt once /3 < oo. Away 
from the symmetry point we have the following: 

Theorem 2.4 If h > h > Q an d jS > ^q{5), then in all translation-invariant, ergodic infinite- 
volume Gibbs states u, 

E ?( ( [S*] 2 ) > 1 - 5, for all r G Z 2 . (2.11) 

Similarly for ) when J2 > ]\ > and /3 > /3o(<5). 

The above results imply the existence of ORLO throughout the part of the quadrant in 
(/i,/2)-plane bounded away from zero. 

2.3. Quantum systems 

While our methods currently seem unable to treat the quantum spin-V2 version of POM that 
was studied numerically in [37], the control of the classical system can be extended to the 
quantum problem with large spin. This is achieved essentially by plugging into the main 
result of [5]. We proceed to introduce the technical aspects of the quantum POM that are 
necessary to state the relevant theorem. 

In the quantum POM, the spins S r are three-component operators S r = (S*,S y r ,Sf) 
satisfying the commutation rules of the Lie algebra su(2), 

[§i,S k r ,] =2i§Xr> (2-12) 

for any cyclic permutation (j,k,£) of (x,y,z). We will work with the (2<S + 1) -dimensional 
irreducible representation of su(2), where S £ {0, \, 1, \, . . .}. This representation is best 
discussed by means of the spin raising and lowering operators 

Sf :=S x r ±iS y r . (2.13) 

The Hilbert space is the linear span of vectors 

\...M r ...) := (g)|M r ), M r = -S,-S + l,...,S -\,S, (2.14) 

r 

on which the operators S z r , S+, S~ act as follows 
S z r \ ...M r ...) = M r \...M r ...), 

S+\...M r ...) = ^S(S + 1) - M r (M r + 1) I . . . Mr + 1 . . .), (2.15) 
S~\...M r ...) = ^S(S + 1) - M r (M r - 1) I . . . M r - 1 . . .). 

t A word of caution: Strictly speaking, the rigorous statement does not rule out the existence of additional 
extremal translation-invariant Gibbs states apart from those above. However, an overwhelming majority of 
the spins in a typical configuration in any such state will be close either to ±ei or to ±e2- 



Long-range order in the plaquette orbital model 9 
The Hamiltonian ii^ on the torus Tjv is then the operator 

K N := -hS~ 2 £ S*S*, - / 2 <S" 2 £ S*S r z ,. (2.16) 

(»"/)* {r,r')z 

The y-component of the spin does not enter the interaction. The scaling by S~ 2 ensures that 
the Hamiltonian is, for each N, bounded uniformly in S > V2. 

The thermodynamical equilibrium is described by means of a linear functional ( - }n,j6 
on the algebra An of all bounded operators generated (via the spectral theorem) by the 
operators S 2 , S+, S~, r G Tjv. Explicitly, 

(A) Nl p= \ m A Aeyi N . (2.17) 

Tr(e-^N) 

Our main result for the quantum system is now as follows: 

Theorem 2.5 Suppose ]\ = Jj > 0. Then for each 5 > 0, there are constants c = c(S) > and 
&o = (5o(S) < °° such that for each 6 > &o(S) and S > cB 2 , we have 

S-\m 2 m 2 ) N ^>l-S, a = x,z, (2.18) 

while 



<i, (2.19) 



and 



< (2.20) 



uniformly in r,r' G Tjy provided N is sufficiently large. 

Note that the bound (2.20) yields the same bound on ( [S z ] 2 [§*,] 2 )n,p because the identity 
I (Ad)nr\ = I (BA}nr\ holds for all self-adjoint operators A, B. It is standard that validity of 
such bounds implies nonanalyticity of the free energy in the appropriate conjugate variables. 
In our case, this will be the function 

f(h) := Urn i z logTr^exp|-/3H N + ^- 2 £ ([S x r ] 2 - [S z r ] 2 ) } J , (2.21) 

where the "external field" h couples to the natural order parameter [S*] 2 — [S z ] 2 . Based on 
this fact, we extract the corresponding result for the asymmetric situations as well: 
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Theorem 2.6 Suppose ]\ > > and let c = c{5) > and jSo = j8o(<5) < 00 be as above. Then 
for each 5, each ft > (5q(5) and S > eft 2 , we have 

S-\[S^] 2 )^>l-25, (2.22) 

while 

s-\m 2 m 2 ) N/fi <s> (2-23) 

uniformly in r, r' G Tjv for N is sufficiently large. A similar result holds for ]j_ > ]\ > with 
indices x and z interchanged. 

We remark that the reason for assuming S > cfi 2 is that, in the underlying proof, we 
use thermal fluctuations to dominate the quantum fluctuations (namely, the effects resulting 
from the non-commutative nature of the relevant variables). 



3. Ground states 



The goal of this section is to prove our characterization of the grounds states on Tjy- As a 
starting point we note the following rewrite of the energy function: 

Lemma 3.1 Let "Hjy denote the torus Hamiltonian. Then for all S, 



n N {s) = 




Proof. Every vertex in has two x-edges and two z-edges coming out of it. Opening up 
the squares in the first two sums, the "diagonal" terms there are easily checked to exactly 
cancel the terms in the third sum. □ 

Proof of Theorem 2.1. First, we claim that 

mmH N (S) = - max{/i,/ 2 }|T N |. (3.2) 

The inequality < is seen by taking S r = e\ or S r = ei, depending on whether ]\ > or not, 
and evaluating Kn(S) for this choice. For the opposite bound we use ]\, > to drop the 
first two terms in (3.1) and conclude 

minH N (S) > -|T N | max {h[S x } 2 + } 2 [S Z } 2 ) > - max{/ l7 / 2 }|T N |. (3.3) 

s SeO(2) 

This also shows that the minimum is attained only by configurations for which 

Sy = S*, for the endpoints r,r' of all a-bonds (r,r') a , oc = x,z, (3.4) 
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and 

7i[S?] 2 + / 2 [S?] 2 = max{/iJ 2 } forallreT N (3.5) 

hold true. 

Let now S be a ground state of Hn- The above findings guarantee that S* is constant on 
any x-plaquette and S z is constant on any z-plaquette. Applications of q> r to the plaquettes 
where the corresponding value is negative defines a configuration S where Sy > for 
a. = x, z and all r £ Tjy. (This is the desired modification of S by applications of the maps 
(^ r ).) Let now r be a lower-left vertex of a z-plaquette and let r' := r — e\ be its neighbor 
to the left. Since S is also a ground state and the edge (r, r') is an x-edge, we have S* = §*,. 
However, the components S z , are both non-negative and since one component of a spin 
determines the other up to a sign, we also have S z = S z . It follows that S r = S r i. Proceeding 
similarly for all pairs of neighbors in we conclude that S r = e for some unit vector 
e e {v e R 2 : v • e { > 0,i = 1,2} and all r e T N . 

It remains to determine the set of vectors e that are admissible at given values of the 
parameters of the model. We have 

n N (~S) = - (h{e ■ df + } 2 (e ■ e 2 ) 2 ) \T N \. (3.6) 

Thus, when ]\ > ]i, we must have e = e\ while e = ei when J 2 > }\. At the symmetry 
point, ]\ = hi an y e W1U gi ve t ne same value. We have thus shown that S is a modification 
of a constant configuration (namely S) of the desired type in all cases of interest. □ 

4. Technical ingredients 

In this section we assemble the technical ingredients needed for the proof of the main 
theorems concerning the phase transition in the model of interest. The proofs will come 
in Section 5. 

4.1. Chessboard estimates 

The proof of the positive-temperature part of the results will be based on the technique 
of chessboard estimates, based on reflection positivity, whose origins go to the seminal work 
of Dyson, Frohlich, Israel, Lieb, Simon and Spencer from the late 1970s. This technique, 
along with a related infrared-bound technology, has proved extremely useful in establishing 
symmetry-breaking phase transitions in various classical and quantum systems with a 
continuous symmetry [18, 14, 16, 17], order-disorder transitions in the Potts and related 
models [29, 10, 3], low-temperature ordering in liquid-crystal models [21, 40]. More 
recently, this technique has also been used to prove phase transitions in systems with 
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highly degenerate ground states without an underlying symmetry [4, 1] including gradient 
fields with a non-convex interaction [7]. The theoretical foundations of this technique are 
well developed already in the original papers [16, 17]; the more recent developments are 
summarized in the lecture notes [2]. 

Consider the model with the Hamiltonian on the torus Tjy with N even. Consider an 
even integer B that divides N and let Ag := {0, 1, . . . , B — 1} x {0, 1, . . . , B — 1} denote 
the block of B x B vertices with the lower-left corner at the origin. Consider a regular 
partitioning of into pairwise disjoint translates of Ag by vectors from {Br:r G T N /g}. 
For r G T N / B , let d r denote the translation by Br. On the configuration space, 

[MS)] r , = S r ,_ Br , r' G T N/ r G T N/B . (4.1) 

We call an event A a B-block event if A depends only on {S r :r G Ag} and we use 3g to 
denote the collection of all B-block events. 

For each A G 3^ we now define a family (0 r (*4): r G T N /g} of translations-reflections as 
follows. First, let A\ denote the reflection of A through the horizontal mid-line { (^j^, y)'.y G 
R/ (NR) } halving the box Ag. Similarly, we use Aj_ to denote the reflection of A through the 
vertical mid-line { (x, ^p): x G R/ (NR) } and .A12 to denote the reflection of A through both 
lines (the two reflections commute and so the order in which they are taken is immaterial). 
Note that A\, Az, An G Ug. For r = rj) G T N /g, we then set 



9r(A) 



r d r 1 (A), for r\j2 even, 

$~ x {A-\ ), for r-i odd and r? even, 

_/ (4.2) 

^ r (^2)/ f° r r i even and rj_ odd, 

k t^jT 1 (^12), for r!,r 2 odd. 



Here := {# r (S): S G ^}. Notice that d r {A) thus depends only on the part of the 

spin configuration in the block Br + Ag. We then have: 

Lemma 4.1 (Chessboard estimate) Suppose > and B > 0. Then for any events 

A\, . . . , A m G $b and any distinct r\,...,r m G T N /g, 



m^n^iAt)) < film J (1 ^(A)) 

V i=l J '=1 L V r€T N/B / 



1 (B/N) 2 



(4.3) 



The punchline of this result is that the probability of a simultaneous occurrence of several 
(often undesirable) events on the torus is bounded by the product of the probabilities of 
events where the individual B-block events have been disseminated — using the maps 6 r — 
throughout the entire torus. The latter quantities are often rather explicitly computable. 

Proof of Lemma 4.1. As already alluded to, the key input for the chessboard estimates is 
reflection positivity of the interaction. We will now define the necessary concept and check 
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the validity of this property. Let P be a plane bisecting a horizontal or vertical line of edges 
of Tjv, i.e., P is of the form either 

{(n + V 2 ,y):y G R/(NR)} U {(« + N/ 2 + i/ 2/ y) :y e R/(NR)}, (4.4) 

with n = 0, 1, . . . , f - 1, or 

{(x,n + V 2 ):x G R/(NR)} U {(x,n + N/ 2 + i/ 2 ) : x e R/(NR)}, (4.5) 

with n = 0, 1, . . . , ^ — 1. The plane has two components and it splits the torus into a left 
half and the right half T^. Abusing the notation slightly, let 6p denote the map on the 
configuration space representing the reflection <R> T^. 

A sufficient condition for the interaction to be reflection positive is that for each such a 
plane P, there are functions g and h = (hi) depending only on {S r : r £ T^} so that 

-H N = g + 6 P (g)+ J^hMhi). (4.6) 

i 

We will now demonstrate that is indeed of this form. Let 

g(S):=-h E S*S$-J 2 E SrS Z r < (4-7) 

(r/)x (r/)z 



and note that g depends only on the spins in T». The collection of functions h will be 
parametrized by the vertices r G that have an edge to a vertex r' € T^. We will use P + 
to denote the set of such vertices r. We set 



Jh S~, if (r, r') is an x-edee, 
h r := { V ll r X ' 6 r G P , (4.8) 

/ J2 S^, if (r, /) is a z-edge, 



where r' stands for the reflection of r through plane P. Then h r 6p(h r ) = JaSyS^,, with a 
depending on the type of the edge (r, r') . A moment's thought then shows that 

-H N = g+9 P (g) + E h r e P (hr) (4.9) 

r€P+ 

and so the interaction is of the desired form. As the a priori measure on the spins has a 
product structure, standard theory (cf [16, Theorem 4.1] or [2, Theorem 5.8]) readily implies 
the desired claim. □ 



4.2. Gaussian calculations 

Through the use of chessboard estimates, the proof of the phase transition will be reduced 
to some tedious but explicit computations of multivariable Gaussian integrals. Informally, 
these can be understood as calculations of spin-wave free energies corresponding to the 
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spin system at hand. The goal of this section is to carry out these calculations and derive the 
necessary estimates between actual partition functions and their Gaussian approximations. 
Throughout we will assume that 

h = k=:J (4.10) 

with / > 0. 

For a given unit vector e{9) := (cos 9, sin 9) G R 2 and A > 0, we define the quantity 

Z N (6,A) := e-^ N \ f e -^ N (s) jj l {lSr _ e{e)l<A} ft v(dS r ). (4.11) 



reT 



reT N 



This is the partition function restricted to configurations within A of a constant configuration 
pointing in direction of the unit vector e(8). Next, we define a function 9 i— > F(9) as follows. 
For each k := (k\, kz) G [—71, n] 2 , let us introduce the quantities 



a± := 1 ± e 



Consider the matrix 



/ 



M(k,6) :-- 



1 
2 



+ \b. 

* 



b± := 1 ± e 



-ik 2 



and p '■= — cos (20). 



(4.12) 



pa_a + 
pb*_b + 




|fl+| 2 + |&_ 



pb-b 


-l 2 +l& 



pb-b* 



+ \ 





\ pb*_b + pa*_a + 

We will see in a moment that det M(k, 9) > and so we may define 

dk 



pa-a^ 

2 



(4.13) 



\b+\ z ) 



HO) -5log( W) + i( " logdetM( M ). 



(4.14) 



The key facts about the quantity F(0) and its relation to Z^(9,A) are the subject of the 
following two claims: 



Proposition 4.2 For any k := (k\,k2) G [— 7r, 7r] 2 and any 6> G [0,27r] we have 
sin 2 (20)sin 2 (fci)sin 2 (fc 2 ) < detM(fc,0) < 16sin 2 (20). 



(4.15) 



Jn particular, F(9) is finite for all 9 with sin(26>) 7^ 0. The function 9 1— >■ F(0) zs periodic with 
period n li, symmetric and continuous on the interval (0, n li), increasing on (0, n l£) and decreasing 
on ( n /4, n li). The infimum ofF is —00 and it is achieved exactly at 9 G {0, n l% tz, 3n /2}- See Fig. 2. 



Proposition 4.3 Suppose that ]\ — h ='■ / > 0. For each r > there are numbers 5 > and 
No < 00 such that if 

/3/A 2 >l PJA 3 <5, (4.16) 




Figure 2. The plot of function 8 \— > F(8) for 8 ranging from to 2n. The key 
fact is that F tends to minus infinity as 9 tends to multiples of n li. However, 
note that F is relevant for the approximation of N~ 2 log 2^(8, A) to within r 
only when 9 avoids the region where | sin(20) | < T — i.e., exactly the vicinity 
of its singularity points. 



then 



< T (4.17) 

holds true whenever N >N$ and | sin(2#) | > t. 

We will first focus on Proposition 4.3 because its proof explains the reasons underlying the 
introduction of the quantity F. The proof consists of a sequence of approximations followed 
by a standard diagonalization of a multivariate Gaussian integral. Let us write each S r as 

S r = (cos(6 + d r ), sin(0 + t? r )) . (4.18) 

On the event entering the integral (4.11), d r = 0(A), so if A is small (which we may assume 
since (4.16) forces A < 5 2 ), then d := (d r ) may be regarded as deviation variables. The 
conversion to a Gaussian integral is performed as follows: 

Lemma 4.4 Suppose A < 1 and define the function 

Gn,oW " <V) 2 sin 2 (0) + \ E (*r~ <V) 2 cos 2 (0). (4.19) 

(r,r')x ( r ' r ')z 

There exists a constant c € (0, oo) such that for every S that is related to # via (4.18), with 
d r e (—71, 7t], and obeys \ S r — e(0)\ < A at allr £ Tjv, 

\n N (S) + J\T N \ - JG N>e (#)\ < c/|T N |A 3 . (4.20) 



Long-range order in the plaquette orbital model 



16 



Proof. With the restriction d r £ (— n, n] and A < 1, the correspondence between S and # is 
one-to-one. The restriction \S r — e(9)\ < A implies d r = 0(A) uniformly in r £ Tln- The 
claim now follows by writing 7^ in the form (3.1), and noting that 

S x r - S x r , = (d r - fy) sin 6 + 0(A 2 ) (4.21) 

and similarly for S z r — S*,. □ 

If we set A' := 2arcsin( A /2) for the maximal angle between S r and e(6) allowed by the 
constraints, we thus have 

where \/27r appears due to our normalization convention for v and where 

Xa>,n{#) ■= E[ !(-A',A')(^). (4.23) 
reT N 

In order to evaluate the expression in (4.22), we will notice the following relation between 
the function Gn,o and the matrices M(k,6). Consider the Fourier-reciprocal torus, 

Tn := { %(ni, n 2 ): n lt n 2 = 0, . . . , N - l}, (4.24) 

and let (T^)+ be its first quadrant, i.e., the collection of k = {ki,k 2 ) £ satisfying 
< ki,k 2 < n. Let M(6) k k / be the N 2 x N 2 matrix indexed by k,k' £ that is block 
diagonal and whose entry is zero unless k,k' £ {fc, + ne\,k + ne 2 ,k + ne\ + ne 2 ) for 
some £ (T^)+, and whose entries for these four values are collected (in the given order) 
in the 4 x 4-matrix M(k,9) defined above. 

Lemma 4.5 We have 

G N ,e(#)= E M(e) r yd r d r/f (4.25) 

r,r'eT N 



where 



M(<0^:=J-r E M{e) kik ,e^- k >. (4.26) 
Proof. We may write G^e as 

= E E 7v(*r-*f+e) 2 , (4-27) 
r6T N e=ei,e 2 

where J rrCl := sin 2 (6>) for vertices r with even first coordinate and := cos 2 (0) for 
vertices r with odd first coordinate, and same for } r /f , 2 and the second coordinate of r. 
Invoking the (discrete) Fourier representation # r = |Tjv| -1 ^ 2 EfceT^ $k e ~ lk r we now get 

Gnjb{*)= E E Wt,(l-e- ike )(l-^ k ' e )J k , k i(e), (4.28) 

k,k'eT N e=e lr e 2 
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where 



kk>(e) : 



Jr,e 

reT N 



e i(fc'-fc)r_ 



(4.29) 



Now, since the couplings ] r>e are 2-periodic in direction e and translation invariant in the 
complementary direction, jfcfc'( e ) will be non-zero only when k! = k or k! = k + ne. In 
these two cases we get 



\ 

Jk,k(e) = 2 and Jk,k+ne(e) 



-cos(20). 



(4.30) 



Consider a fixed vector k e and consider the restriction of the sum in (4.28) to 

k, k' £ {£, + ttci, k + 71^2/ ^ + tt^i + tt^}- Let a±, b± denote the quantities in (4.12) for 
this k. The quadratic form corresponding to e := e\ is then described by the matrix 



/ 



\ a - \ 
pa*_a+ 







pa-a + 
\a+\ 2 













pa-a* + 

2 



\ 



(4.31) 



pa_a + \a+\ ) 



while the contribution corresponding to e := is described by the matrix 



/ I&-I 2 

|fc_| 2 
pb*_b+ 

V pb*_b+ 



pb-b% \ 
pb-b 



\K\ 
o 







'+\ 



J 



Adding these contributions together, we obtain (4.25)-(4.26) with (4.13). 



(4.32) 



□ 



Proof of Proposition 4.3. Let 1(8, A) denote the integral in (4.22). Our goal is to evaluate 
1(9, A) to within multiplicative correction of the order e°^ T ^ TN L If it were not for the 
indicator Xa',n(^)> the integral would be Gaussian; unfortunately, as Gjv,e(#) depends only 
on the differences d r — d r i, it would also diverge. We will therefore have to treat the indicator 
with some extra care by deriving suitable upper and lower bounds. 

First, for any A > 0, 

Xa',nW <exp{^5Jc^ 2 |T N |-^6/ £ tf 2 }, 



(4.33) 



reT t 



where c' := sup| A | <2 (AVA) 2 £ [1, oo). Substituting this bound into 1(6, A), we now scale d r 
by -y/j6j and use Lemma 4.5 to diagonalize the quadratic form in the exponent to get 



1(0, A) < (^-^iTNle^/c'A 2 !^! Yl 



ke(A) + y/det{A + M{k,e))' 



(4.34) 
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The A in the denominator now regularizes the contribution of the k = mode. We thus get 
the bound 

lQg ^ A) < C/3/A3 + ^/3/c'A 2 - F N (6,A), (4.35) 

where 



F N (0,A) := UogtfJ) + E logdet[A + M(M)]. (4.36) 



N\ 



Since F^(6, A) is, in the limit N — > oo, finite and continuous in and also obviously larger 
than F(6) for any A > 0, it suffices to choose A so that the first two terms on the right 
of (4.35) are strictly less than t. In particular, taking A := A, it suffices to choose 8 so 
that (c + \c')8 < t and pick N so that F N (6,A) > F(6) for N > N and all d, to get 
N~ 2 log Zn(6,A) + F(9) < t whenever BJA 3 < 5. 

The requisite lower bound is derived by a change-of-measure argument. Consider the 
Gaussian measure underlying the upper bound above: 

P A (d^):=el TN l^^ A )-^ G ^Wexp{--/3/ £ tf 2 } f\ -^L (4.37) 

2 reT M reT N V 2n 



Using E A to denote the corresponding expectation, for any A > we clearly have 

Z N (6,A) > e- c ^ A3 \ T ^-\ T ^ e ^E A (xA',N)- (4-38) 

Since F^(6,X) will tend to F(8) when N — >• oo and A J, 0, uniformly in N with | sin(0) | > T, it 
will suffice to show that, for any A > small enough, Ea(^a' / n) 1 ^' Tn ' is near one uniformly 
as N — > oo. To this end, we first use chessboard estimates — exactly in the setting described 
above with B := 1 — for the Gaussian measure Fa to derive 

Ea(Aa',n) = E A ( n ^I^KA'}) > E[ F a(I^| < A ; ). (4.39) 

reT N reT N 

Next the Chebyshev inequality yields 

P A (|*,| > A') < P A (|tf,| > A) < ^1 = (4.40) 

where we used that A ; > A. The variance of a Gaussian variable increases when we make 
the corresponding quadratic form in the exponent smaller (as a matrix). As this quadratic 
form is bounded below by A/3 / times identity, the thus get 

P A (|tf r -0| <A) >l-^i- (4-41) 

To derive the matching lower bound from (4.38), we thus need to choose A small so that 
Fjq(6, A) is already close to F(8) but such that A/3/A 2 is still large. A specific choice that will 
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work is as follows: Pick 5 > so small that, for some No, we have \F N (d, V$) - F(0)| < \ for 
all N > Nq and all 6 with | sin(26>) | > T. In addition, assume that also c5 — log(l — y/S) < 5. 
The choice A := y~8 and the bounds (4.38) and (4.41) then yield the desired lower bound on 
log Zn(6, A) once both conditions (4.16) are satisfied. □ 

Proof of Proposition 4.2. Assume that k is such that a±, b± ^ and recall the definition 
of p. First, notice that detM(fc, 8) = when p = ±1. The case p = 1 is easily checked by 
multiplying the matrix M(k, 6) by the vector (1, —1, —1, 1); the case p = — 1 is checked using 
the vector (1, 1, 1, 1). Since det M(k, 8) is an even quartic polynomial in p that is divisible by 
(1 — p 2 ), it can be written in the form 

det M(k,0) = (1 - P 2 ){A - p 2 C). (4.42) 



An explicit computation yields 
1 

16 



A := ^(|fl-| 2 + \b-\ 2 ) (K| 2 + |fc_| 2 ) (|fl_| 2 + |fc+| 2 ) (K| 2 + \b+\ 2 ) (4.43) 



and 



C:= ^(|fl_| 2 |fl+| 2 -|&_| 2 |fc + | 2 ) 2 . (4.44) 



Obviously, under the conditions a±,b± ^ we have A > and C > 0. Moreover, since 
| fl ±|/ |b±| < 2, we have A < 16. Therefore 

|fl_| 2 |a + | 2 |&-| 2 |& + | 2 < A _ c < A _ p 2 c < A < 16 _ (4 45) 

The left-hand side equals sin 2 (fci) sin 2 (/c2) and so (4.15) follows by plugging the above 
inequalities into (4.42). 

To get the second part of the claim, let g(x) := (1 — x)(A — Cx) be the function appearing 
on the right-hand side of (4.42). Theng'(O) = —A - C < and g'(l) = -A + C < 0. But 
g is quadratic and so it is strictly decreasing throughout [0, 1). In particular, g(l) < g(x) for 
all x G [0, 1). Using this in (4.14), the desired claims follow. □ 



4.3. Good and bad events 

For the proof of our key estimates, we will follow, as in the previous sections, the general 
scheme developed in [4, 1] and further discussed in [2, Sect. 6.4]. For a positive constant n 
and an even integer B that divides N, define the good B-block events Q x and Q z as follows. 
First set 

Q° x := n {S: \S r -e x \ <rj}, (4.46) 

reA B 
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and similarly for with e\ replaced by ej_. Now let Ag denote the even-sub lattice vertices 
in Ag; these are the lower-left corners of either x or z-plaquettes. We use A B to denote the 
set of those vertices in Ag that are lower-left corners of x-plaquettes, A| = Ag \ Ag. For any 
set A C Ag, let <p A (S) be the composition of q> r for all r € A. As these maps commute, the 
order of the composition is immaterial. With the help of these notations, we now set 

&== U 9a@1), a = 1,2. (4.47) 

AcA° 

These are the good B-block events; the requisite fcfli event is defined by 

B:= (g x ug z ) c . (4.48) 
For a B-block event A define the quantity 

(B/N) 2 



pp(A) := lim 



V fd., ,„ / 



(4.49) 



where N is taken to infinity along the even multiples of B. The limit exists by standard 
subadditivity arguments. Note that this quantity is a limiting version of the objects on the 
right-hand side of the chessboard estimates. A key input for our proofs is the observation 
that pp(B) is small: 

Proposition 4.6 For each n > and e > there exists 6q < oo such that for any & > Bq there 
exists B for which the bad event B defined using these n and B satisfies 

Pp(B) < e. (4.50) 

In order to show that the bad event is unlikely to occur, we will need to further partition 
it into several subevents. First, consider a number A > — to be evetually chosen in 
dependence of B so that the conditions (4.16) hold true — and let us isolate the configurations 
where the corresponding component of two neighbouring spins somewhere in Ag differ by 
more than a quantity proportional to A: 

Be:= U U { SeS: l S ?" S r'l >T^}- (4-51) 

" =X ' Z (r,r% 
r,r'eA B 

These will be the configurations with too much energy — hence the subscript "E". The 
remaining "bad" configurations will be collected in the set 

Bsw:=B\B E , (4.52) 

where "SW" designates the method — a spin-wave calculation — that will be used to estimate 
the contributions to this event. The following claim provides a key structural information 
on the configurations contained in £>sw : 
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Lemma 4.7 Suppose that n and A satisfy the inequalities 

A < n 2 < 1. (4.53) 
Then for each S £ £>sw/ £/iere is a unit vector v £ R d and a set A c A B swc/z £/iflf 

|(f> A (S)) f -t>| < |, r£A B . (4.54) 

Proof. Fix S £ B sw and define A C A° g by taking AHAJ = {r £ A B : < 0}, oc = 1,2. 
Consider the configuration S := q> A (S). Clearly, S* > for all r £ A|, a = x,z. The proof of 
the main claim now comes in two steps. First we will show that both components of all S r , 
r £ Ag, are at least v 7 ^ in absolute value, where 

f ==£-!• W-55) 

This will be used to ensure that > everywhere. Then we will rerun part of the argument 
to show that S is, actually, to within 0(A) of a constant configuration. 

Let £ be as above and abbreviate c := (16B) -1 . Notice that the bounds (4.53) imply 

C + j < \ cA < ^ and > 2cA. (4.56) 
First, we claim that 

S £ S sw => [S?] 2 > £, r £ Ag, a = x,z. (4.57) 

Focusing our attention on on a = 1, assume that [S x ] 2 < £ at some r £ Ag and derive a 
contradiction with S £ £>sw- Notice that S £ £>sw implies that 

|S?-S£|<cA (4.58) 

for any r, / £ Ag connected by an a-edge. If r' is a neighbor of r over a x-edge, then 

[S/] 2 < iVI + cA} 2 < £ + 2cA, (4.59) 

where we used that £ < \ and cA < \. On the other hand, if r' is connected to r by a z-edge, 
we can proceed via the second components to get again 

[S x r ,] 2 = 1 - [S z r ,] 2 (4.60) 

< 1- [S z r ] 2 + 2\S z r ,-S z r \ (4.61) 

< £ + 2cA. (4.62) 

Examining all pairs of nearest neighbors along a shortest path from r to any r' £ Ag and 
using that this path has at most 2B edges, we conclude 

[S x r ,] 2 < C + (2B)2cA = Z, + j, r'e A B . (4.63) 
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Notice that to maintain the argument from (4.59) for the proof of the inequality [S*„] 2 < 
[S*,] 2 + 2c A for every step (r 1 , r") of the iteration, we actually need to invoke that £ + - < \. 



At the same time, since y 1 — (C + §) > 5 > 2c A, the triangle inequality implies > on 
each z-plaquette, and thus everywhere in Ag. In addition, this allows us to compute further: 



\~S r > -e 2 \ 2 = [S x r ,] 2 + (l - y/l - [S^]2j (4.64) 

< [S^] 2 + [S^] 4 < 2[S*,] 2 <2^+j)=2C+^ = V \ (4-65) 

where we used that 1 — \A — x < x for < x < 1. This would imply that S G and thus 
S G Q z , in contradiction with the assumption that S G £>sw C 

Having proven the bound (4.57), we can use the fact that y/T, > 2cA to imply that S* > 
for all r G Ag. Combining now the positivity of components of Sy with (4.57), we can 
improve the bound on the difference of the 1-components of spins over a z-edge (r,r') z as 
follows: Suppose without loss of generality that S*, > S z and write 

(4.66) 



Pr' 5 r 



[s z r ,V-Ji-[s 



< < (4-67) 

1 - [S5] 2 VW 



In conjunction with | S z , — S z \ < cA, this implies 

|S/ - S r \ < cA f^ ^ 4cA ( 4 - 68 ) 

where we used that 1 — £ > i. As (2B)4cA < -, the desired claim now follows with v := Sq 
and A as above by examining a path of minimal length between the origin and any other 
site in Ag. □ 

The previous lemma allows us to further partition £> sw as follows. Let v\,...,v n denote 
the unit vectors representing the complex n-th roots of unity, V£ := eii 1 ^), £ = 1, . . . ,n, 
where n : = |_^-J + 1 is the smallest integer for which nA > 2n. Defining 



B (i) 



(J <p A ({S e B sw :Vr e A B \S r -v e \ < A}), (4.69) 



ACA U B 

it follows that 

#sw = U B®. (4.70) 

Having in mind that A <ti n and that, by definition, £>sw H (Qx U Qz) = / w e notice that 
= whenever the distance of vg from the points ±ei, ±e2, is less than n — A. The set 
function A 1— >• pg(*4) is subadditive (see [2, Lemma 5.9]) and so we have 



p 6 (£)<p 6 (£ E ) + £p 6 (^). (4.71) 

l=\ 
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It remains to derive suitable estimates on pg(£>E) and pg(£>g^). An indispensable ingredient 
will be the following lower bound on the full partition function: 

Lemma 4.8 Fix x > such that sin(2r) > x and let 5 and No be as in Proposition 4.3. If A and B 
satisfy the conditions (4.16) and N > No, then 

{Z N ,p) 1/lTNl >eV- F V-\ (4.72) 

where F is the free energy introduced in (4.14). 

Proof. By restricting the integrals over the spins to the set where \S r — e(x) \ < A, we get 

Z N #>e^ T ^Z N {x,A). (4.73) 
From here the claim follows by invoking the bound (4.17). □ 



Lemma 4.9 Suppose that n G (0, n /4) and, given x > with sin(2r) > x, let 5 be as in 

■TV 

4 



Proposition 4.3 and assume that A and (5 satisfy the conditions (4.16) and that A < n 2 < 1 



Then 

Pp{B E ) < (2B) 2 e c i B2 - c ^ A2/B2 , (4.74) 
where c\ := \ log(27r) + F(t) + x and cj_ := 1/512, and 

P/3(4w) < 2 y *e-^-*)-HT)-^\ t = 1/ ... /n . (4.75) 

Proof. Consider the disseminated event B E := f]r6T N/B ^j-(-^e)- We need to derive an upper 
bound on the constrained partition function Zjy r(£?e), which is given by the same integral 
as the full partition function but only over configurations in the event Sg. Since on £>e, each 
configuration has an "energetically charged" edge in each block in {Ag + Bt: t G T N /s}, 
Lemma 3.1 tells us 

1 A 2 

n N (S) > -(— ) J(N/B) 2 -J\T N \, SeB^. (4.76) 

To account for the entropy, we note that there are altogether 2B(£> — 1) positions where the 
"energetically charged" edge can occur in each translate of Ag. Lemma 4.8 now gives 

^(B3< B/N > 2 =(^gi)) <B/N)2 (4.77, 
< ( 27r )8 2 /2 (2B) 2 e - W AV512+[FM + TlB^ (4 . 78 ) 

where the factor (2tt) b2 ^ 2 comes from performing (unconstrained) integrals of the spins. 
The limit N — >■ oo along multiples of B now yields (4.74). 
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Next we consider the event B^ := nreT N/B $r(£>sw)- We only need to focus on £ 
satisfying \^£ — ^fk\ > f] — A for = 0,1,2,3 because otherwise the event is void. Here 

we notice that, for each S G ^sw' there is a unique A C UteT N/B (Bt + Ag) such that the 
configuration S := <P/^(S) is within A of the vector vg. The number of distinct A associated 
with a single S is 2^ N ^ B ^ 2 ^ 4: . As both the a priori measure and the Hamiltonian are invariant 
under q> A , this yields 



Z N;/3 (S^ ) < 2( N / B ) 2 / 4 Z N {Pf, A) . (4.79) 

Hence, 

m^(^) (B/N)2 = (4-80) 

<2 1/4 exp{[-F(^) + F(t) + t]B 2 }. (4.81) 

Now we invoke Proposition 4.2 and apply the monotonicity and periodicity properties of F 
to infer that F(£^) > F(n - A) for any £ for which ^ 0. (Th is is where we need 
< A < n < n /i.) The limit N — > oo (along multiples of B) finishes the claim. □ 



Proof of Proposition 4.6. Fix e > and 77 £ (0, rc/4) and set x to be a positive number such that 
sin(2r) > x and 

F(n/ 2 ) >F(t)+t. (4.82) 

This is possible because 1— > F(0) tends to minus infinity as 6 I 0, see Fig. 2. Let (5 be related 
to this T as in Proposition 4.3. We will henceforth link A and B to /3 > via 

A: =j 6-n an d B:=2Llog j 6j. (4.83) 

Notice that this choice of A will make the bounds in Proposition 4.3 true once B is sufficiently 
large. Let Bq be now a value such that for all B > Bq these bounds hold, the condition (4.53) 
in Lemma 4.7 is satisfied, and the inequalities 



V - A > I (4.84) 



and 



(2B) 2 e c i B2 - c ^ A2 / B2 + ^L 2 1/ *e-^- F ^-^ B2 < e (4.85) 

hold true. This is possible by (4.82) and our choices (4.83) — and the fact that the quantities C\ 
and C2 depend only on x. 

We claim that (4.85) implies the desired bound pp(B) < e. Indeed, using that n < 47 7a 
for sufficiently small A, it suffices to insert into (4.71) the bound on pp(B^) and the uniform 
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bound on pp(B^) from (4.74-4.75). As a result we get the estimate of pp(B) by the left-hand 
side of (4.85). (We also used that F is strictly increasing on (0, n l^) and applied (4.84).) By 
(4.85), pp(B) is thus less than e for any 6 in excess of Bq defined above. □ 



5. Proofs of main results 

5.1. Classical model 

We are now ready to prove our main results for the classical system. The key inputs are 
chessboard estimates alongside with the bound in Proposition 4.6. 

Proof of Theorem 2.3. The standard line of reasoning leading to the proof of phase coexistence 
in the present context is based on the observation that distinct types of good events are 
unlikely to occur in the same configuration. A formal statement is as follows: 

Assume that B > 4 and n < 2 sin(7r/ 8). For any e > there exist 8 > and No such 
that if pp(B) < 5 and N > N , then ft Ni p{6 f (g x ) n r '(£ z )) < e for all r,r' e T N/B . 

The proof of this fact is an application of a Peierls-type contour argument. Fix r,r' G TT N / B , 
r 7^ r' , and let y(r) denote the set of all A C T N /b such that both A and T N / B \ A are 
connected and A 3 r. We use 9A to denote the set of all vertices outside, but adjacent to A. 
For A £ y(r), consider the event 

Ra--={E ls°0 r >i|8A|} 7 (5-1) 

r'edA 

where 6 r (B) = ft r (B) with 6 r denoting the shift by Br. Recalling also the notation cf\, resp., 
£72 for the shift by 2e\, resp., 2ei, we now claim 

6r{G x )ne r ,(g z )cz [j (r a u (J (o-j(R A )uo--\R A )) r ). (5.2) 

Aey{r) 7=1,2 

To prove this, consider a configuration S G 6 r (G x ) Pi6 r r(Q z ) and let A' denote the set of all 
s £ T N / B for which there is a nearest-neighbor path that starts at r, ends at s, and visits only 
blocks Bs' + Ag, where 6 s i (Q x ) occurs. Clearly, r' A' so we may use A to denote the set of 
all s £ T N /g such that every nearest neighbor path from s to r' contains at least one vertex 
of A' — i.e. A is A' with all of its "holes" filled. 

Our task is to show that 

S e R A U (J (aj(R A ) U o--\R A )) r . (5.3) 

7=1/2 
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If S G Ra, then we are done, so let us suppose that S ^ Ra- Under this condition more 
than 4/5 of all blocks corresponding to vertices in 9 A are good but, since they are not part 
of A, they are of type Q z . Then, however, more than 2/5 of Q z -blocks in 9 A are adjacent to 
A in one coordinate direction — say it is the first one — and thus 1/5 of them are adjacent 
in either positive or negative coordinate direction. Observing that S £ 9 r {Q x ) H 6 r+e .(Q z ) 
implies that cr 1 " 1 (S) G r (B) — and similarly for the opposite direction — we get (5.3). 

With (5.3) on hand, we now perform a standard version of the Peierls argument combined 
with chessboard estimates. A key input is the fact that, for some c G (1, oo), 

|{A G y(r) : |9A| = n}\ < c n . (5.4) 

Using inclusion-exclusion and cy-invariance of ji^^, it now suffices to show that for some 
c < oo and No < oo, 

jMRa) < cl 3A ^(£)l aA l /5 , A G y(r). (5.5) 

To this end, let No be such that, for all N > No, the probability on the right hand side of 
(4.49) for A := B is at most [2pp(B)\ N . Now, cover Ra by the union over all subsets of 
3 A with |3A|/5 elements where the bad event occurs. The number of such partitions is at 
most 2 9A ; the probability of each occurrence is estimated by [2pp(B)] |3A|/5 . This proves 
(5.5) with c := 2 6 ^ 5 and thus the above claim. 

Applying this claim alongside the fact that, on a good event, the spins are oriented along 
one of the coordinate directions proves (2.8) and (2.9). To get also (2.10), one just follows 
literally the argument proving the main result in [6]. □ 

Proof of Theorem 2.4. This is a consequence of Theorem 2.3 and a convexity argument. 
Let £ x = £q denote the Wenzel-Janke plaquette energy (1.5) for the type-x plaquette at 
the origin. Let <3(}i,}2) denote the set of all translation-invariant, ergodic Gibbs states. 
The convexity of }\ i— > logZ^ then implies that once ]\ < }[, for any ]i G S5(/i,/2) and 
}i' G <&{][, h) we have 

E^(£ x ) < ]Ey(£*). (5.6) 

Now, at ]\ — J2 —'• J and /3 > /3o, Theorem 2.3 guarantees the existence of a ^ G (S (/,/), 
such that 

E^(S X ) >A-S. (5.7) 

(Here S may differ from the one used in the statement of Theorem 2.3.) Thus, for all ]\ > ]i 
and y. G <5(/i, we have 



(5.8) 
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As £ x is the sum of four terms of the type S*S*„ which are less than one, we must have 

E H (S?S$) >1-S (5.9) 

for all nearest-neighbour pairs (r,r') x of type x. This forces E f ,(|S*|) > 1 — 5 and thus 
E ?i ([S?] 2 ) > 1 - 25. As 5 is arbitrary, this proves the claim. □ 

We will also give the formal proof of absence of magnetic order: 

Proof of Theorem 2.2. It is easy to check that every Gibbs measure is invariant under the action 
of any <p r As each spin belongs to one plaquette of type-x and one plaquette of type-z, it 
follows that the distribution of Sy is symmetric around zero, for any r and any oc = x, z. This 
implies the claim. □ 



5.2. Quantum model 

Our set of results for the quantum model will be derived by an application of the general 
theory developed in [5] whose main conclusion can be found in Theorem 3.7 of [5]. 
This theorem says roughly the following: Whenever a quantum-spin model satisfies the 
conditions of (quantum) reflection positivity, and the classical system admits a proof of 
phase coexistence at a positive temperature by means of chessboard estimates, then the same 
phase coexistence occurs in the quantum system provided the magnitude of the quantum 
spin is sufficiently larger than the inverse temperature squared. 

First let us check that the prerequisite concerning the quantum reflection positivity is 
satisfied. Since we are using reflections in planes bisecting edges of T^, this is proved by 
the same argument as in the classical case, except that we need to write all operators in a 
basis in which their matrix elements are all real valued. This is satisfied automatically in the 
representation (2.15) in which S Z Y , Sf are real and so is S* = ^(S+ + S~). 

All we have to do is thus adapt the proof for the classical model to plug into Theorem 3.7 
of [5]. We begin by introducing the formalism of coherent states that the whole connection 
is based on. Consider the space 

C 25+1 that 

carries the corresponding 2<S + 1-dimensional 
representation of su(2). Let Q be a vector on the unit sphere §2 in R 3 that is described by the 
spherical angles 6 and (p. Then we set 

|0> : = E ( S + m) /2 [ cos ^)] 5+M [sin(^/ 2 )] 5 - M e^-^ |M>. (5.10) 

Abusing the notation slightly, whenever O := (O r ) is a collection of such vectors, we will 
denote the corresponding product state by |Q) := (g) r |fi r ). 

The coherent states have a number of remarkable properties of which relevant for us are 
particularly those listed in Sect. 2.1 of [5]. Here we will only need the notions of the lower and 
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upper symbols. Given a linear operator A on <S> r eA C 25+1 , the lower symbol is the function 
H i-> (A) n on (S 2 ) A such that 

(A) n := (Q|A|Q). (5.11) 

The upper symbol is, in turn, a function O i— » [ A] q such that 

2<S + 1 X 



A =\^T) y (S2)A dn[A] n |n)(Q|, 



where dO := Eire A dn r with dO r denoting the uniform measure on §2 with total mass 4zr. 
The upper symbol is not necessarily unique, so we use the notation [A]q to denote any 
version thereof. The upper and lower symbols are two natural classical approximations of 
the quantum Hamiltonian, so we need to check that they are reasonably close: 

Lemma 5.1 Consider the operator %^ in (2.16). Then there is a constant c = c(Ji, J2) < 00 and a 
version of the upper symbol [Hn]q such that for each G £ (§2) Tn / 

|<#N>n-[#N]n| <||T N |. (5.13) 

Proof. Let 'Hn(S) denote the classical Hamiltonian corresponding to (2.16) — obtained by 
replacing the operators S r by vectors S r G §2 and dropping the normalization by S 2 . Then 
the desired bound will hold once we verify that 

(%)n = H N (H) and [ft N ] n = (1 + 1/<S) 2 H N (Q) (5.14) 

for some version of [Hn]ci- This is in turn shown by noting that is multilinear in the 
operators S — which means that for each term in the sum the corresponding symbols take 
the form of a product — and by the fact that (S r )o = SCl r and that [§ r ]n : = («5 + l)O r is a 
version of the upper symbol for S. The bound (5.13) now follows by the fact that 7^(0) is 
bounded by a constant times |T]v|, uniformly in O and S > V2. □ 

Our next step will be to verify the conditions required by [5, Theorem 3.7] for the 
corresponding classical model. The problem here is that the link through the coherent 
states naturally leads to f/zree-component classical spins. Consider the events Q x , Q z and 
B that are defined as follows. For a three-component spin configuration S = (S r ) with 
S r = (Sy, S^, Sy), let S xz denote its projection onto the xz-plane scaled to have a unit length. 
For a B-block event A C (S 2 ) Tn , let 

A := {S: S xz e A} n p| {S: \S y r \ < A} (5.15) 

reA B 

denote its natural extension to three-component spin configurations. This immediately 
defines the events Q x , Q z , and B. However, B does not cover the complement Q x U Q z ; for 
that we will also need the event 

B y := |J {S:\S l l\ > A}. (5.16) 

reA B 
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Then, obviously, B := {G x ^Qz) c satisfies B = B U B y . In order to plug into our calculations 
for the two-component spin, we will need to prove the two lemmas. 

Lemma 5.2 Suppose A < V2. Then there exists a constant C3 < 00 such that, for any B-block 
event A such that AC\B E = 0, its counterpart A defined as above satisfies 

Proof Consider the expectation of the disseminated event A in the torus measure fix a for 
the three-component model. We will write this expectation as Z N) p(A)/Z N/ p. Pick S G A 
and define S := S xz . First, by the very definition S belongs to the disseminated event A. 
So in order to compare the expectation for the three-component spins with that for the two- 
component spins, we will need to control the change in the Hamiltonian and the a priori 
measure under the map S 1— > S. 

To treat both Z^^{A) and Zjv,^ in a unified fashion, assume that \S v r \ < Sior some 5 < A 
and all r e T N . Let S := S xz . Note that 

S a r = S&r = S a r + 0{5 2 ), a = x,z. (5.18) 

From the fact that A n B E = it follows that 

(§ a r - S*,) 2 = (S« - S a r ,) 2 + 0(S 2 A) (5.19) 
and so we have 

n N (~s) = n N (s) + j £ [s y ] 2 + o(s 2 a)}\t n \. (5.20) 

reT N 



The a priori measures are related by v(dS r ) — y 1 — [S y ] 2 v(dS r )dS y . 

We will now derive bounds on the partition functions Zjy a (A) and Zjv,^ in terms of their 
two-component spin counterparts Z^^{A) and Z^r. For Z^ r ^(A) we set 5 := A and note 
that Hn(S) > %n(S) + 0(A 3 )/3/|T;v|. Integrating out the components S y gives 

Z N ^(A)<Z N 4A)(ne°^f N ^ (5.21) 

For a lower bound on Z N /J g, we first note that Z N ^ > Z N ^(B^ n Q'), where Q' z is the event 

G' y = fl {S:\S y \<5} (5.22) 

reA B 

with £ := A 2 . Applying (5.20), the Hamiltonians now differ by 0(A 3 ) /|T;v| and the integral 
over S y , r G Tjv, now yields a term A 2 per site. This shows 

Z^>(A 2 e°( A3 )W) |TNl Z^(^). (5.23) 
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Combining the upper and lower bounds and applying the subadditivity bound p«(Bg) > 
1 — p«( Be) we get the desired claim. □ 

Lemma 5.3 There are constants c±, C5 £ (0, 00) such that 

~ / *\ B 2 r2 p -c 5 /3/A 2 

max{p^(H E ) /P ^(H y )} < (c 4 A- 2 e^ A3 ) j _ (5e) (5-24) 

Proof (Sketch). We will instead estimate the objects P/3(Be \ -By) and pg(By \ Be); from these 
the claim will follow by invoking subadditivity A 1— >■ pp(A). Consider the partition function 
Zjv,^(Be \ By) and let S be a configuration contributing to the sum. As the second component 
of all spins are small, the partition function is bounded as in the two-component spin case. 
The result is 

Zn^We \ By) < e-WI T ~l [ 2y /to] I Tn > [B 2 e-^ 2 ] {N/B) \ (5.25) 

Combining this with the lower bound on Zjq a from the previous proof now gives the desired 
bound for P^(Be \ By). 

Concerning the event By \ Be, as this is disjoint from Be, we just note that, in this case, 
that we can bound 

H N (~S) >n N (S)+JA 2 (^y + 0(A 3 )J\T N \ (5.26) 

and then proceed as in the previous lemma. The desired suppression now comes from the 
second term on the right hand side of the last display. □ 

We are now ready to prove our main result on quantum systems with large spins: 

Proof of Theorem 2.5. Key to the formalism of [5] is the operator Q.4 associated with the 
event A as follows: 

Qa-= — / dO |0)(Q|. (5.27) 



4tt ) Ja 

Notice that if A is a B-block event, then Q.4 behaves as identity on the part of the Hilbert 
space outside Ag. The map A i— >• Q.4 is countably additive; in particular, if A\, . . . , A n form 
a partition of the probability space, then Qa x + • • • + Qa„ = 1- 

Consider now the operators Qg , Qg and Qs- Notice that these operators are invariant 
under reflection of the box through any mid-plane (because so are the events they arise 
from). Let d t denote the shift on Tjv by Bt, where t £ Tjy/g. Let Z, := ^, where c is as in 
Lemma 5.1. Since B < C\y/S and p^(B)e^ +C2 ^^ v ^ is small — of course, for & large, thanks 
to Lemmas 5.2 and 5.3, Proposition 4.6 and the choices (4.83) — for some absolute constant 
C\, C2 £ (0, 00), Theorem 3.7 of [5] tells us that for some e > small 

(QsW<e (5-28) 
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and 

(Q eti ^(l-Q^ z )> N ^<e. (5.29) 
for any t Xl tj_ G T N/B . 

In order to process these to the bounds in the statement of the theorem, let us note that, 
if O is a bounded operator generated by with oc — x, y, z, and r G (tB + Ag) U (t'B + Ag) 
and [0]a is its upper index for which we set 

7 :=sup|[0] n | (5.30) 
n 

and 

7 a := sup||[0] n | : O £ OtQa UOfGaj, oc = x,z, (5.31) 
then (5.29) and (5.31) imply 

(0) Nt p < 3ej + max{ 7x ,7 z }. (5.32) 
Indeed, writing 

(0) N f = ^— |TWl / dO [Ofo (0|e-^|0), (5.33) 

we can estimate [0]n by y a on 6tG a U # f /(? a , and by 7 on the complement of these. The 
positivity of (Q|e _ ^ N |n) permits us to convert the resulting integrals to expectations of 
the kind (5.28-5.29). 

To see how this applies in a specific situation, consider O := S~ 4 [Sr] 2 [S^,] 2 for r,r' 
belonging to distinct translates of Ag. We then get 

[0] n = [n^f^f + OOS- 1 ), (5.34) 

which is 0(n + <S _1 ) on U a=z , z (6>f£? a U 6 t /Q a ) and 0(1) otherwise. As a consequence, for 
some c < 00 and N is sufficiently large, 



s- 



$f [^] 2 >n,a ^fo + S^ + e). (5-35) 



The other cases needed to establish (2.19-2.20) are checked analogously. Once (2.19-2.20) 
are proved, the bound (2.18) follows by the symmetries of the model and the fact that 

<s~ 2 Ejs?] 2 = 1 + 0(1/5). □ 

Proof of Theorem 2.6. Fix / > and throughout this proof let denote the Hamiltonian for 
Ji = J2 — }■ Consider the operator 

£%:=S~ 2 £ S x r S x r „ cc = x,y,z (5.36) 

<r/>« 
r/eT N 
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and notice that increasing ]\ above the common value / amounts to adding the term 
(/ ~~ Ji)&n to Hn- Key to the proof is to show that 

1 Tr(e~^ N+hS N) 
f(h):= lim — log^ s '- (5.37) 



satisfies 

d 



dh 



-f(h)\ h=0 > 2(1-3). (5.38) 



Indeed, h i— >• f(h) is convex (by Holder inequality) and so (5.38) implies the same bound for 
all h > 0. This in turn shows that, for N sufficiently large and r with both coordinate even, 

S- 2 (£ r ) Nrfi >2(l-6) (5.39) 

where £ r is quantum counterpart of the plaquette energy (1.5). From here the desired claims 
follow along the same argument as in the classical case. 

We thus have to show (5.38). To this end consider the torus events 

:= |q: ^ 0*0*, > (2-c5)|T N |], a = x,z. (5.40) 

<f/>« 

A straightforward application of chessboard estimates shows that, for each 5 > 0, there is 
j S 1 = fatf) such that for > fa and S > cfi 2 , we have 

1 

where — > as N — >■ oo. Next we note that, by Jensen's inequality 

Tr ( e -^ N +^) > ^ 2 ^+^ ' TnI ^ dQ e -/ 3 ^N>n+M2-^)|T N | / (5 42) 

where we already applied that (o^}ci > (2 — <5)|Tn| for O G *4 a - Theorem 3.1 of [5] now 
shows that 

e -/3(#N>n > (n|e"^ N |n)e" c ^l TN l /v ^ (5.43) 
for some c < oo. Substituting this into (5.42) yields 

Setting a := x and applying (5.41), the bound (5.38) follows once 
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6. Concluding remarks 

6.1. Neel vs orientational order 

As already mentioned, Wenzel and Janke determined in their numerical experiments that 
the model exhibits a Neel ordering of the plaquette energies; see Figs. 2(a-b) of [37]. 
Explicitly, energy was found low on the z-plaquettes and high on the x-plaquettes in a 
sample obtained by multiple updates of a configuration using Glauber dynamics with the 
Metropolis rule. We wish to point out that this ordering is a direct consequence of the ORLO 
established rigorously in the present work. Indeed, if the spins are with high probability 
aligned with ei — i.e., pointing north or south — then the z-plaquettes will have energy 
E r ~ —4 while the x-plaquettes will have energy E r ~ 0. In addition, the mixed plaquettes 
will settle at energy E r ~ —2. This is consistent with Fig. 2(d) of [37]. 

The state with the spins aligned with e\ will have the roles of the x and z-plaquettes 
interchanged, giving the energy distribution again a Neel type order resemblance. 
Notwithstanding, the physical significance of a Neel order is unclear given the period-2 
nature of the interaction. And, in fact, matters seem to look quite different when instead 
of £ r we work with more natural quantity, 

8 r (S) := [S? - S? +e f + [S? - S? +e f 

r s a 3 _ s a 3 1 2 , [c«4 _ c«4 1 2 / fi v 

^[ D r+e 1 D r+e 1 +e 2 i ^ V°r+e 2 D r+e 1 +e 2 i ' v - 1 ; 

which is the plaquette energy normalized to vanish in all ground states. Although we do not 
see any reason why a strong Neel ordering should be exhibited by these plaquette energies, 
it would be perhaps of some interest to rerun the numerical experiments at higher precision 
to check this fact numerically. 



6.2. Enhanced sampling 

The samples of actual configurations shown in Fig.2(d) of [37] deserve one more comment. A 
cursory look at the figure reveals some level of orientational order in the e2 _s pi n direction — 
which is consistent with our mathematical results — but a more careful analysis uncovers an 
apparent statistical discrepancy. Indeed, most of the spins point down in the figure although 
the plaquette-flip symmetries (q> r ) of the Hamiltonian, which can be applied independently 
at all even locations, indicate that about half of the z-plaquettes should be pointing up! It 
is thus somewhat surprising that in Fig.2(d) of [37] only 7 such plaquettes out of the total 
of 25 have upward-pointing spins; the others are clearly pointing down. As the number of 
upward-pointing plaquettes is well approximated by a binomial distribution, the probability 
that this happens is 

( f ) ^5 * 0.014, (6.2) 
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i.e., the configuration in the figure will typically appear only once in about 70 samples! 
We take this as a possible indication that, despite judicious methods of simulations, the 
configuration may not have fully equilibrated at the time the snapshot was taken. 

Turning this observation into a positive statement, one can try to use the plaquette-spin 
flips to accelerate the convergence of the computer sample to equilibrium. Notice that, even 
at moderate temperatures, any single-spin update rule will have considerable difficulties to 
overcome the energy barrier associated with changing the orientation of an entire plaquette. 
The dynamics would naturally mix faster if an occasional flip of an entire plaquette — 
by an application of one of the maps (q> r ) — is incorporated into the stochastic dynamics. 
This would result in an algorithm reminiscent of the Swendsen-Wang method for sampling 
configurations in the ^-state Potts model [35]. It could be expected that this enhancement 
would result in a substantially better performance of the simulations. An algorithm of this 
sort has been recently attempted in the context of the orbital compass model [38], although 
there the cluster flips have to be performed along entire lines of sites which makes them 
very non-local. 

6.3. Correlation decay 

Our mathematical argument establishes rigorously long-range order in the system. 
Nonetheless, we do so without giving any bound on the decay of (truncated) correlations. 
Since our argument is based on contour methods and suppression of long-wavelength 
part of the spin-wave decomposition, we tend to believe that the correlations generally 
decay exponentially fast in any translation-invariant, ergodic Gibbs measure for this system. 
However, we have not been able to find a rigorous argument in the vain. Again, it would be 
of interest to see if this question could be addressed by numerical methods. 
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